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Artifactual log-periodicity in finite size data : 
Relevance for earthquake aftershocks 

Y. Huang^'^, A. Johansen^, M. W. Lee^, H. Saleur^ and D. Sornette^'"^'^ 

Abstract. The recently proposed discrete scale invariance and its associated 
log-periodicity are an elaboration of the concept of scale invariance in which 
the system is scale invariant only under powers of specific values of the 
magnification factor. We report on the discovery of a novel mechanism 
for such log-periodicity relying solely on the manipulation of data. This 
"synthetic" scenario for log-periodicity relies on two steps: (1) the fact that 
approximately logarithmic sampling in time corresponds to uniform sampling 
in the logarithm of time; and (2) a low-pass-filtering step, as occurs in 
constructing cumulative functions, in maximum likelihood estimations, and 
in de-trending, reddens the noise and, in a finite sample, creates a maximum 
in the spectrum leading to a most probable frequency in the logarithm of 
time. We explore in detail this mechanism and present extensive numerical 
simulations. We use this insight to analyze the 27 best aftershock sequences 
studied by Kisslinger and Jones [1991] to search for traces of genuine log- 
periodic corrections to Omori's law, which states that the earthquake rate 
decays approximately as the inverse of the time since the last main shock. 
The observed log-periodicity is shown to almost entirely result from the 
"synthetic scenario" owing to the data analysis. From a statistical point 
of view, resolving the issue of the possible existence of log-periodicity in 
aftershocks will be very difficult as Omori's law describes a point process with 
a uniform sampling in the logarithm of the time. By construction, strong 
log-periodic fluctuations are thus created by this logarithmic sampling. 



1. Introduction 

Log-periodic corrections to scaling have been one of 
the most exciting topics in the physics of rupture in 
recent years [Saleur and Sornette, 1996; Saleur et ai, 
1996a, b]. If present, these corrections might make the 
prediction of major events a less hopeless task [Ani- 
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frani et ai, 1995] than previously thought [Geller et at, 
1997]. Log-periodicity is a characteristic behavior of 
systems exhibiting discrete scale invariance (DSI) [Sor- 
nette, 1998], either geometrical or spontaneously gener- 
ated by the dynamics. Like all symmetries, DSI could 
provide a powerful new concept in the study of cracks, 
rupture, and disordered systems in general. 

DSI is an elaboration of the concept of scale invari- 
ance, which has been broadly applied in the geosciences 
[Scholz and Mandelbrot, 1989; Barton and Pointe, 1995; 
Turcotte, 1997; Sornette, 2000]. The word "fractal" was 
coined by B.B. Mandelbrot to describe sets which con- 
sist of parts similar to the whole and which can be de- 
scribed by a fractional dimension. This generalization 
of the notion of a dimension from integers to real num- 
bers reflects the conceptual jump from translation in- 
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variance to continuous scale invariance. More recently, 

the generalization of the notion of dimension, accord- 
ing to which the dimensions or exponents are taken from 
the set of complex numbers, has been shown to describe 
the interesting and rich phenomenology of systems ex- 
hibiting discrete scale invariance, associated with log- 
periodic corrections to scaling [Saleur and Sornette, 
1996; Saleur et ai, 1996a, b; Sornette, 1998]. Discrete 
scale invariance is a weaker kind of scale invariance in 
which the system or the observable obeys scale invari- 
ance only for specific choices of magnifications, which 
form, in general, a countably infinite set of values, all 
powers of a fundamental scale factor A. This property 
can also be seen to include the concept of lacunarity of 
the fractal structure. 

Log-periodic oscillations have been reported in acous- 
tic emissions prior to rupture [Anifrani et al, 1995; 
Sahimi and Arbabi, 1996; Johansen and Sornette, 1998; 
Sornette, 1999], to rock bursts in deep mines [Ouil- 
lon and Sornette, 2000], in precursory seismic activity 
before large earthquakes [Sornette and Sammis, 1995], 
and in anomalous ionic concentration prior to the Kobe 
earthquake [Johansen et al., 1996, 2000a]. The status 
of log-periodicity is the clearest for rupture of hetero- 
geneous materials for which a sub-harmonic cascade of 
MuUins-Sckerka instabilities has been recognized as the 
probable underlying mechanism [Sornette et al., 1996; 
Huang et al., 1997]. The application to earthquakes is 
controversial [Gross and Rundle, 1998] as only a few 
case studies have been reported. In this context, a ma- 
jor difficulty consists of the selection of the appropriate 
space-time window to construct timc-to-failure func- 
tions [Bowman et al., 1998; Brehm and Braile, 1999]. 
Another problem stems from the scarcity of seismic data 
relevant to this problem since only large (and thus rare) 
earthquakes exhibit precursory seismic activity observ- 
able solely on relatively large foreshocks [Reasenherg, 
1999]. Notwithstanding these limitations, further in- 
vestigation of the possible existence of log-periodicity 
is important for its potential use (1) in constraining 
the mechanisms controlling earthquake triggering, (2) 
in focusing possible forecasting skills based on the mea- 
sured accelerated seismicity, and (3) in providing in- 
direct measurements of the stress drop ratio and the 
seismic efficiency [Lee, 1999]. 

It is fair to say, however, that the existence of log- 
periodic corrections in rupture is not proven beyond 
doubt: They have typically been observed in systems 
with very noisy data, and so far, we are still lacking a 
nontrivial analytically solvable model where they could 
be shown to exist rigorously. By nontrivial, we mean 



a realistic model of rupture without preexisting hierar- 
chy which could be treated analytically to prove log- 
periodicity. In that respect, the situation is not unlike 
critical phenomena, where the possibility of power law 
singularities was hotly debated before Onsager [1944] 's 
solution of the Ising model. This must be contrasted 
with the preexisting hierarchical geometry, such as the 
hierarchical diamond lattices or tent-like hierarchies or 
Cantor sets used by Newman et al. [1995] and Saleur 
and Sornette [1996]; Saleur et al. [1996a, b] for which 
the log-pcriodicity derives rather trivially by measuring 
the fractal dimensions of the underlying geometrical ob- 
ject [Sornette, 1998]. 

It is therefore of vital importance to have a greater 
understanding of the relation between apparent log- 
periodic corrections and noise. After the initial period 
of enthusiasm, the suspicion has grown, in particular, 
that the consideration of cumulative quantities, like the 
Benioff strain in foreshock sequences, might well intro- 
duce spurious log-periodicity in the data, despite the 
lack of any physical DSI. The purpose of this paper is to 
examine several scenarios where apparent log-periodic 
corrections are generated by the interplay of noise and 
the method of data analysis. We report on a basic 
artifactual source generating log-periodicity, which is 
a combination of a sampling that is periodic in a log 
timescale and some type of low-pass filtering such as in 
constructing cumulative functions, maximum likelihood 
estimators, and de-trending. 

In the analysis of power law data, sampling is more 
often regular (sometimes even periodic) on a log scale 
than might be initially expected. For instance, in our 
previous analysis of the Kobe earthquake the chemi- 
cal concentrations reported by Tsunogai and Wakita 
[1995, 1996] were sampled in a nonuniform way because 
the measurements were performed after the earth- 
quake. Specifically, the chemical composition of ground- 
water prior to the earthquake was analyzed by collecting 
bottles of mineral water from wells close to the epicen- 
ter for which the production dates could be identified. 
Since the depletion process of bottles in stores can be 
approximated by a Poisson process with some rate fi, 
this implies that the number of bottles with a produc- 
tion date t prior to the date of the earthquake tc is pro- 
portional to exp[— /X {tc — t)], quantifying an accelerated 
sampling in the time tc — t to the earthquake. This is 
not exactly a sampling that is uniform in the variable 
In {tc — t); it will nevertheless be approximately so in 
finite time intervals and, as a consequence, may lead to 
the effects discussed in this paper. 

Another example occurs naturally when analyzing 
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aftershocks. Omori's law for aftershocks [Utsu et ai, 
1995] states that the rate Q{t) of earthquakes foUowing 
a niain earthquake decays as t"^, where t is the time 
from the main shock and p is close to 1 : 

Q{t)ocl/tP. (1) 

It is equivalent to having the probability of an after- 
shock in the time interval [t, t+dt] proportional to dt/f^ . 
Alternatively, the occurrence time t for an aftershock is 
a random variable with a power law probability distri- 
bution. When analyzing a given empirical aftershock 
time series, the problem of measuring the rate precisely 
becomes crucial. To compute the rate, one has to choose 
a binning interval. A constant interval is not satisfac- 
tory since Omori's power law implies that events cluster 
densely at the beginning of the time series and become 
very sparse at long times. A binning interval varying 
with time is the most natural, but it will, when opti- 
mally chosen, realize once again a periodic sampling in 
log scale. We will discuss this in great detail in Section 
3. 

The spurious log-periodicity generated by such sam- 
pling will usually become amplified when considerat- 
ing cumulative quantities, which act as a low-pass fil- 
ter. This is quite unfortunate, since taking cumulative 
quantities is a very natural procedure in many cases. 
For instance, taking the cumulative of the empirical 
rate of aftershocks appears to be a natural well-defined 
approach, which has the advantage of being nonpara- 
metric and amounts simply to counting the number of 
events up to a running time t. 

The paper is organized as follows. In Section 2 we 
analyze in detail the effect of cumulating quantities on 
log-periodicity. In Section 2.1.3 we show that when the 
binning of a power law is logarithmic in time, the cumu- 
lative distribution creates splendid and totally spurious 
log-periodicity. Section 2.2 applies this insight to real 
aftershock sequences: Using this cumulative method, 
we find that it is not possible to distinguish any real 
log-periodicity from the synthetic mechanism in the 27 
best aftershock sequences studied by Kisslinger and 
Jones [1991] and in the Loma Prieta, Landers, and 
Northridge aftershock sequences. Section 2.3 shows 
that, conversely, cumulating can destroy a preexisting 
log-periodicity, when the binning is performed regularly 
and nonlogarithmically. Section 3 turns to another an- 
alyzing method using maximum likelihood. We again 
show that the synthetic mechanism for log-periodicity 
operating in the cumulative method is present here, and 
we quantitatively demonstrate its main properties. Ap- 
plied to this set of aftershock sequences, we find that 




Figure 1. Mapping between a uniformly distributed 
random variable x and a random variable t with power 
law distribution: (a) linear scale for t and (b) log scale 
for t. 

most of the observed log-periodic signal can be ex- 
plained by the synthetic mechanism. However, there 
seems to remain a residual signature, possibly stemming 
from a physically based log-periodicity, which can not 
be accounted for by the synthetic mechanism. This ev- 
idence is not sufficiently strong to conclude definitively 
with high statistical significance. Section 4 summarizes 
our findings. 

2. Cumulative Distributions and Their 
Fluctuations 

Omori's law (1) and, more generally, any continuous 
law can be obtained from a one-to-one mapping to a 
uniformly distributed random variable x (without loss 
of generality, we take x € [0,1] when not specified oth- 
erwise) : 

x^ [ Q{t')dt' , (2) 
Jo 

which is the cumulative distribution function (cdf) ob- 
tained from the condition of "conservation of proba- 
bility" P{x)dx = Q{t)dt under a change of variable. 
Consider the special case p = 1 . This leads to 

P{x)dx ^ dx ^ —dt ~ c?lni , (3) 

which highlights the mapping between the linear scale 
in X and the logarithmic scale in t. 

The mapping (2) is shown graphically in Figure 1 . A 
given value xo on the x axis corresponds to a point to on 
the t axis, through the transformation (2) represented 
by the continuous curve in Figure la. If xo is a random 
variable uniformly distributed on the x axis, to will be 
a random variable with a power law distribution. A se- 
quence of equally spaced x is mapped onto a sequence of 
t with increasing separation following a power law p(t) . 
Note that we are not sampling a power law function, 



4 



HUANG ET AL. 



rather the power law probabihty density function (pdf) 
results from the spacing between the times t of the af- 
tershocks. When t is plotted on a log scale (Figure lb), 
we see that if there is a cluster of points near xq, there 
will be a cluster of points near Iq. When the clustering 
of points in x is periodic on the x axis, the clustering is 
periodic on the \nt axis; thus periodicity in the cdf of 
X is mapped to log-periodicity in the cdf of t. 

2.1. Theory 

2.1.1. Shot-noise spectrum for cumulative dis- 
tributions of finite data sets. To understand the 

effect of constructing cumulative quantities, we con- 
sider the cdf for the variable x uniformly distributed 
in [0,u]. The cdf is by definition cdf(a;) = x/u defined 
for < a; < M. Here wc add the interval size variable 
u to investigate how the fluctuations, in particular the 
most probable ones, may depend on the size of the sam- 
pled interval. Let us assume that a finite number N of 
values X are sampled in the interval [0,u], and let us 
construct the empirical cdf cdfe(a;) from this data set 
by simply counting the number of values less or equal to 
X and normalizing by A'^. We expect to get an approxi- 
mate straight line with step- like fluctuations around it. 
These "shot noise" fluctuations are described by 



1 ^ 

Cdfe(x) = ^ X] ~ 



(4) 



where the Heaviside function: H{y) = for y < and is 
1 otherwise. To study the fluctuations around the theo- 
retical cdf(a;) = x/u, we estimate the Fourier transform 
of the residual R = cdfe(a;) — x/u, which reads 



R = -l- 

iu> 
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(5) 



The flrst term on the right-hand side gives the usual 
fluctuation spectrum = N~^ui~'^ of "red noise." 
The other term in the r.h.s. is only present for finite in- 
terval size u. In the presence of a finite u, the spectrum 
presents correction factors with amplitude propor- 
tional to (^e"^'" — l)/iu)u and its modulus squared. These 
corrections describe the possible existence of oscillations 
with angular frequency u) decorating the linear depen- 
dence. This shows \is that the finite size u of the interval 
introduced a nonzero most probable angular frequency. 
Here, by most probable, we mean the frequency cor- 
responding to the largest peak of the spectrum. It is 



defined as the angular frequency that maximizes 

2 



1 



= — 2 [1 - COs{iOu)] , 



(6) 



which is the solution of 

tan(ciJM,/2) = Luu/2 . 

Thus we find uju/2 « 'Sn/2. Since ui 
the most probable frequency 



(7) 

27r/, this gives 
(8) 



which is inversely proportional to the length u of the 
interval. This constitutes our first prediction. 

Prom the structure of (5), we see that the ampli- 
tude of the leading correction to asymptotic spectrum 
N~^LU~^ is proportional to that is, we expect the 
most probable oscillation decorating the linear cdf x/u 
to have a typical amplitude proportional to 1/ \/]V. This 
constitutes our second prediction. 

2.1.2. Spectrum of cumulative variables that 
are solution of a stochastic ordinary differential 
equations. To show the generality of these results, 
we can generalize the problem within the framework of 
stochastic ordinary differential equations and study the 
"cumulative" variable X defined by 



dX/du = aX + ri{u) 



(9) 



where 77 is a white noise defined by {rj{u)ri{u')) = S{u — 
u'), where the angle brackets mean that an ensemble 
average is taken. Eq. (9) allows for a deviation from 
the simple linear cumulative law by introducing a first- 
order autoregressive correlation (the aX term). The 
solution of (9) is 



X{u) 



Jo 



t]{t) exp[a(u — r)]dr 



(10) 



Eq. (10) shows that X is the cumulative of r], weighted 

by the exponential term exp[a(M — r)] (equal to 1 for 
a = 0). The white noise can be written as 

v{t) = jduj e--^(a;) , (11) 

which defines f]{Lo)- The stochastic variable is char- 
acterized by the following covariancc: {f]{uj)ff {u')) = 
5{u> — uj'), where the asterisk denotes the complex con- 
jugate. Putting (11) into (10), wc get 



X{u) = j du fu{co)e 



,i4>(ui) 



(12) 



LOG-PERIODICITY IN EARTHQUAKE AFTERSHOCKS 



5 



where 



tu — a 



(13) 



The most probable angular frequency is again obtained 
by maximizing |/„(a;)p. Let us take the limit of small 
a for which the noise dominates (in the power law map- 
ping, this corresponds to an exponent p close to 1). 
Then (13) simplifies into 



1 



ICO 



|/„(u;)|2 = (2/c^2)[1-cos(wm)] 



(14) 



(15) 



This retrieves the same solution (8) as above for the 
shot noise problem. Note that the l/w^ spectrum for a 
signal with an infinite range is modified by the cosine 
correction observed for a finite range. 

These simple calculations can be explained intu- 
itively as follows. Consider an uncorrelatcd Gaussian 
noise. In the frequency domain the spectrum of Gaus- 
sian noise is flat (white), since all frequency components 
arc equally probable. Let us call G(lu) the white spec- 
trum of the measurement error g{x). G{u) does not, 
however, extend to all frequencies and is bounded by the 
measurement range of x. The lowest angular frequency 
<^min is determined by the span of x, corresponding to 
one half cycle covering the whole interval (Nyquist fre- 
quency). The highest angular frequency is determined 
by the distance between successive Xi. 

When noise is integrated, as occurs in constructing a 
cumulative distribution, the resulting noise is smoother 
and its spectrum is reddened with more power in the 
low-frequency band. This results from the fact that 
the spectrum of g{x')dx' is G{uj) [Stratonovich, 
1967]. 

Regardless of the precise shape of G(w) (we allow 
here for departure from white noise), the combination of 
the lower angular frequency cutoff and the increase 
implies that lo~^G{u)) presents a well-defined peak near 
Wmin- The left side of the peak is determined by the na- 
ture of the frequency cutoff, while the right side is con- 
trolled by the w"^ decay resulting from the integration. 
This angular frequency peak corresponds to the most 
probable oscillation in the signal, i.e., the frequency at 
which the spectrum is maximum. In other words, the 
integration of the white noise in finite time series cre- 
ates a low-frequency periodic oscillation decorating the 
average trend of / g{x')dx' . 

We thus obtain the following general result: for 
approximately uniformly distributed random variables 




Figure 2. Mapping of a periodic function in x under 
the transformation (2) for p = 1 (solid line), p = 2 
(circles), and p = 0.01 (stars). 



and finite samples, the most probable oscillation deco- 
rating the cumulative distribution is inversely propor- 
tional to the range u of the variables. In addition, the 
amplitudes of the oscillations decay on average as the 
square root of the sample size. This is therefore a finite 
size effect. 

Slutzky [1937] similarly found that taking sums and 
differences of discrete time series, i.e., filtering-smoothing 
them, introduces oscillations in otherwise completely 
noisy data. See also Zajdenweber [1976] for simulations 
with an application to economy. 

2.1.3. Application to log-periodicity in cumu- 
lative power law distributions. With the change of 
variable u -^Int obtained in (3), forp = 1, (9) describes 
a power law X ~ decorated by noise. The previous 
calculations in the linear variable u translates directly 
to the variable Int, predicting a most probable log- 
periodic component with log frequency 1.5/ ln(i2/ti), 
where [ti,t2] defines the interval in the variable t. 

When the power exponent p is not exactly 1, the 
mapping curve in Figure lb will not be a straight line, 
and then the mapping of a periodic function in x under 
the transformation (2) will not yield an exact periodic 
function, as shown in Figure 2. 

In Figure 2 when the power law exponent p is 1 (solid 
line), the function is periodic in Int. For p = 2 (circles) 
the function is not exactly periodic in In t. Owing to the 
rapid decay of the events cluster closer to t = 
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than for the 1/t decay. In contrast, for p < 1 the events 
cluster at larger t more than for the 1/t decay, which 
leads to the oscillations shown by the stars in Figure 2 
for p = 0.01. For more reasonable values closer to 1, 
we observe significant regular oscillations which are ap- 
proximately log-periodic. Their significance is assessed 
by the existence of a significant peak in the spectrum. 
The approximate log frequency peak will be higher in 
the case of p ^ 1 than in the case for p = 1, but the 
difference is not very strong for the empirically relevant 
cases for which the apparent exponent docs not deviate 
very much from 1. In the case of aftershocks, Omori's 
law is typically characterized by exponents in the range 
[0.5,2.0] [Kisslinger and Jones, 1991]. Similarly, the 
accelerating Benioff strain of precursory earthquake ac- 
tivity prior to a main shock is typically characterized 
by a power law {t,-; — t)^P with an exponent p close to 
0.5 [Bowman et al, 1998]. 

Therefore the generation of log-periodicity in the 
case of logarithmically sampled power laws relies on the 
following ingredients: (1) noise and sampling fluctua- 
tions (shot noise) provide all kinds of frequency com- 
ponents; (2) integration filters out high-frequency com- 
ponents, enhances low-frequency components, and pro- 
duces a spectrum only perturbed for the lowest 
frequencies for the finite range in Int of the data set; 
and (3) there is a direct mapping from periodicity to 
log-periodicity for power laws with exponents close to 
p = 1 since data points are close to being uniformly 
distributed in \nt scale. 

Summarizing, we expect to observe log-periodic os- 
cillations from the integration of noise in power laws : 
(1) the (log-) frequencies of the log-periodic oscillations 
should have a power spectrum given by w"^; (2) the 
most probable frequency should scale as the inverse of 
the range in log scale of the analysis; (3) the most prob- 
able frequency should increase when the power law ex- 
ponent p deviate from 1; and (4) the amplitude of the 
oscillations should scale as the inverse square root of 
the number of data points (central limit theorem). 

In Section 2.2 we present numerical tests of these pre- 
dictions performed on synthetically generated samples. 

2.2. Numerical Simulations 

Our numerical simulations are based on the hypothe- 
sis that the noise is uncorrelated. Empirical realizations 
of the noise in real data sets may be of a different na- 
ture. However, we believe that our results hold by and 
large as our mechanism relies mainly on the existence 
of an increasing power spectrum at low log frequencies 
bounded from below by finite size effects. 
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Figure 3. Circles represent the cumulative distribution 

of a realization of synthetic aftershocks with Omori's 
law {p = 1) obtained with (16). The solid line is a fit 
to the synthetic data using a power law. The verti- 
cal axis is the normalized cumulative number of events. 
The horizontal axis is time t (the same unit as the real 
sequence). 

A sequence of event occurrence times ti for i = 
1,..,N are generated within the fixed time interval 
[to, tu], with a power law rate 1/t^ by the transformation 
method [Press et al., 1992]. In this method, one starts 
from a uniformly distributed random variable x £ [0, 1] 
and generates the corresponding random variables t us- 
ing the transformation (2), i.e., 



t 



i-p 



(16) 



on each of the realizations of the x variable. In this 
procedure we specify the number N of points x in [0, 1], 
generate them in this interval using a random number 
generator, and apply the transformation (16) to each of 
these values to obtain the sequence ti for i = 1, ..,N. 
The cumulative distribution of a typical realization (cir- 
cles) is shown in Figure 3. 

In our analysis, we first determine the best power 
law fit to the synthetic data of Figure 3 and study its 
residue, shown in Figure 4. We then perform a spectral 
analysis of the residue to test for the presence of possi- 
ble periodic oscillations. Since the data points are not 
evenly spaced in the Int scale, we use the Lomb method 
[Press et al., 1992]. Recall that for unevenly sampled 
data the Lomb method is superior to fast Fourier trans- 
form (FFT) methods because it weights data on a "per 
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Figure 4. Oscillations around the power law (solid 
line) of the normalized cumulative number of events 
(circles) in Figure 3. The horizontal axis is Int. 
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Figure 5. The spectrum of the residue shown in Figure 
4. The horizontal axis is frequency co. The vertical axis 
is the normalized Lomb spectral power. 



point" basis instead of a "per time interval" basis. Fur- 
thermore, the significance level of any frequency com- 
ponent can be estimated quite accurately [Press et al., 
1992]. 

In order to determine the presence or absence of a 
periodic signal in the residue, it is necessary to quan- 
tify the statistical significance of the peaks found in 
the spectrum Pjv(w). We use the null hypothesis that 
the noise realizations at distinct times are uncorrelated 
Gaussian distributed, which in the Lomb method leads 



to PAr(a;) with an exponential probability distribution 
with unit mean. Thus, if we scan M independent fre- 
quencies, the probability that none has a spectral peak 
larger than z is (1 — e~^)^ . The false alarm probability 
of the null hypothesis is then 

P(peak >z) = l-{l- e"^)^ . (17) 

A small value for the false alarm level indicates a highly 
significant periodic signal. In general, M depends on 
the number of frequencies sampled, the number of data 
points, and their detailed spacing. Usually, M is of 
the order of a few N, the number of data points (for 
example in the Interactive Data Language (IDL) im- 
plementation of the Lomb method, M = iY). Figure 
5 shows the spectrum of the residue presented in Fig- 
ure 4. One observes a main peak near lo = 1.5, with 
height (with normalized spectral density) around 8.5. 
Since the number of data points here is N = 31, using 
the standard choice M = N = 31 gives the false alarm 
level 0.006 for this peak, indicating that the probabil- 
ity that this peak results from chance is 0.6%, which 
corresponds to a confidence level better than 99%. 

In the following, we will consider only the most signif- 
icant peak of our synthetic samples. We first run 1000 
simulations with N = 100 data points with to = 0.1 and 
tu = 100. Using the estimation of the false alarm proba- 
biHty (17), we find that 99.4% among the 1000 synthetic 
time series have a significance level better than 99% for 
the presence of a periodic component. The distribution 
of these angular frequencies is shown in Figure 6 and is 
compatible with the theoretical prediction ~ given 
by (15). 

We now test the prediction (8) that the most prob- 
able angular frequency is inversely proportional to the 
range ln(ti„ax/imin) in Int of the data set. In our tests, 
we set tiniii = 1 and sweep imax- For each value of fmax 
we construct 100 runs to get the average angular fre- 
quency. Since the angular frequencies are distributed 
according to a power law (approximately w^^), the av- 
erage is not the most probable value. However, both 
quantities have very similar distributions, and there is 
a simple proportionality relation between the average 
and the most probable value. We analyze the average 
angular frequency for computational simplicity. Figure 
7 shows the average angular frequency as a linear func- 
tion of l/[ln(fmax/imin)]- 

As seen in Figure 2, when the power law exponent 
p is not equal to 1, the mapping between periodicity in 
the X and log-periodicity in the t variable is not exact : 
the clustering of data points either at small t {p > 1) 
or at large t {p < 1) effectively shortens the logarith- 
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Figure 6. The cumulative distribution of the most 
significant angular frequencies of the log-periodic oscil- 
lations obtained from the cumulative power law distri- 
bution. The straight line is the theoretical prediction 
of (corresponding to slope -1 in the cumulative dis- 
tribution) and is included to guide the eye. 
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Figure 7. The average angular frequency as a function 

of l/[ln(tmax/tmin)]- 
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Figure 8. The most probable angular frequency w as 
a function of the power law exponent m. 

mic range and thus increases the effective value of the 
most probable frequency. We verify this prediction by 
generating, for each value of p, 20 runs that allow us 
to get an estimation of the log frequency distribution 
from which the most probable log angular frequency is 
extracted. We verify in Figure 8 that the most probable 
log angular frequency increases as p departs from 1. 

We now test the prediction that the amplitude of 
the log-periodic oscillation should be proportional to 
AT"^/^, where A'' is the number of data points. For 
each value of N the amplitude is averaged over 20 runs. 
Figure 9 shows the dependence of the log-periodic os- 
cillation amplitude for N going from 100 to 5995, with 
the best fit ~ N~^''^^ in excellent agreement with the 
prediction. 

Although log-periodic oscillations from integration of 
noise tend to be low log frequency, corresponding to one 
to two oscillation cycles, our simulations show that the 
probability of observing oscillations that span up to 2.5 
periods is still around 10%. 

2.3. Tests on Real Aftershock Data 

Recent numerical simulations and a mean-field anal- 
ysis have been performed on models of earthquake after- 
shocks in one-, two-, and three-dimensional Euclidean 
lattices, using both continuous elasticity and discrete 
cellular automata modeling [Lee, 1999; Lee and Sor- 
nette, 2000]. Owing to the interplay of power law 
stress corrosion and the threshold nature of rupture, 
the average stress is found to decay in a punctuated 
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N 

Figure 9. The scaling relation between the amplitudes 
of the log-periodic oscillations and the number of data 
points N. The slope is -0.49. 

fashion after a main shock, with events occurring at 
characteristic times increasing as a geometrical series 
with a wcll-dcfincd multiplicative factor (a signature of 
log-periodicity) which is a function of the stress cor- 
rosion exponent, the stress drop ratio, and the degree 
of dissipation. These results, found independent of 
the discrete nature of the lattice, raise the possibility 
that log-periodic structures may also be observable in 
real aftershock data [Lee, 1999]. If confirmed, their 
existence would provide constraints on their underly- 
ing mechanism as well as on the stress corrosion ex- 
ponent, the stress drop ratio, and the degree of dissi- 
pation. Another possible mechanism comes from the 
hierarchical structure of the fracturation of the crust 
with the coexistence of an approximate self-similarity 
together with characteristic scales associated with the 
thicknesses of geological and rheological layers [Ouil- 
lon et al, 1995, 1996]. These considerations have mo- 
tivated us to search for log-periodic oscillations around 
Omori's power law of the normalized cumulative num- 
ber of events for 27 aftershock sequences in southern 
California, 1933-1988. The log-periodic components in 
the oscillations are extracted using the Lomb method 
and compared with those from synthetic aftershock se- 
quences generated with the same parameters as the real 
aftershock sequences. 

In order to avoid possible bias in the selection of the 
aftershock time scries we decided to study a sufhciently 
large and high-quality data set that has been studied 
previously by an independent team, namely, the 39 af- 



tershock sequences in southern California covering the 
time from 1933 to 1988, compiled by Kisslinger and 
Jones [1991] using the Reasenherg [1985] algorithm. In 
addition, we studied more recent large aftershock time 
series from the Loma Prieta earthquake (October 1989 
in northern California), the Landers earthquake (June 
1992), and the Northridge earthquake (January 1994). 

Kisslinger and Jones [1991] used these sequences to 
study the properties of aftershocks in southern Cali- 
fornia, focusing on Omori's power law l/V^ and possi- 
ble correlation between the exponent p and geophysi- 
cal properties like surface heat flow. They found that 
the temporal behavior of aftershocks is generally well- 
described by the modified Omori's law. The p values 
range from 0.69 to 1.809, with mean value of 1.109 and 
standard deviation 0.246. Among the sequences, six 
values were graded "poor," 10 "fair," and 24 "good,", 
on the basis of the relative size of the standard error to 
the p value itself with a standard error < 0.1 p being 
good, > 0.2 p being poor, and fair being in between. 
The number of aftershocks in one sequence ranges from 
12 to 1437, and the duration of the sequences ranges 
from 7 days to around 1000 days. Note that the sizes 
of the samples are adjusted by the chosen magnitude 
cutoffs which have been varied to test the robustness of 
our results. 

In order to test for the possible existence of log- 
periodic oscillations decorating Omori's power law we 
need the 1/V power law to be a good approximation 
in the first place. We thus discarded the poor data 
sets (4, 7, 8, 26, 34, 35) as defined by Kisslinger and 
Jones [1991]. We then examined the other sequences 
very carefully before using them, and some more se- 
quences were discarded for some obvious reasons: 

1. A sequence is severely incomplete in the early 
stage as seen from the plot of the temporal distribution 
of the events (4, 7, etc., confirming Kisslinger and Jones' 
fits). 

2. There is a sudden change of properties in the mid- 
dle of the sequence due to some large aftershocks which 
probably caused their own aftershocks (19, 21, 37, 38). 
This is explicitly illustrated for the two sequences 19 
and 37 (for sequence 19, p = 5.8 ±9.0 from a fitting pro- 
gram in a commercial software system Matlab using the 
Levenberg-Marquardt algorithm, while Kisslinger and 
Jones' result is p = 1.113 ± 0.415. For other sequences, 
our results are close to their results. For sequence 37, 
there are two subsequences with very different p, one 
is close to 1.1, the other is close to 2.0 according to 
Kisslinger and Jones [1991].). If the sequence is suffi- 
ciently long, we truncate it and keep the part before the 
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distortion for the analysis. Otherwise, we discard it. 

3. The number of events is too small (say < 15) or 
the duration is too short (say < 1 day), for example, 
events 16 and 18 last 0.21 days, events 13 and 39 last 
15 days. 

Of the 39 aftershock sequences studied by Kisslinger 
and Jones [1991], 12 were thus excluded from our anal- 
ysis. For the other 27 sequences we first change the 
magnitude cutoff above which the aftershocks are kept 
in the analysis. This allows us to determine the best 
magnitude cutoff for which the oscillations decorating 
Omori's power law for the normalized cumulative num- 
ber of events are the most regular. This selection in- 
troduces a slight bias in favor of the presence of log- 
periodicity, which in fact strengthens our negative con- 
clusion obtained by this approach. In any case, similar 
results are obtained by choosing the magnitude thresh- 
old used by Kisslinger and Jones [1991]. The regularity 
and strength of the oscillations are quantified by the ra- 
tio between the sum of the squares of the errors (RSSE) 
in the time domain resulting from the power law fit and 
from a power law with log-periodic oscillations added on 
it. The presence of log-periodic oscillations is also quan- 
tified by the significance level of the main peak of the 
oscillation spectrum obtained from the Lomb method. 

We illustrate this procedure on the 1933 Long Beach 
earthquake sequence, with a main shock magnitude of 
6.3, which contains 571 events with minimum magni- 
tude 2.0 and total duration of 2576 days. For this se- 
quence, using a magnitude cutoff of 4.5 provides the 
strongest signal: we find the smallest correction for the 
origin of time of the aftershock sequence (a fit to the 
modified Omori's law l/(t + Cy gives a very small C, 
an exponent p close to 1 and a high RSSE). The Lomb 
periodogram of the oscillations has also a very signifi- 
cant frequency peak. The cumulative number of events 
is shown in Figure 10, in which the solid line is the fit 
to a pure power law. The residual of this fit is shown in 
Figure 11 and is fitted to a simple cosine function. The 
corresponding spectral analysis is shown in Figure 12. 

Synthetic tests were performed for this magnitude 
cutoff of 4.5 by generating 200 numerical sequences with 
exactly the same parameters as determined by the sim- 
ple power law fit. Figure 13 shows that the real se- 
quence (defined as the intersection of the two lines) ex- 
hibits a log-periodic frequency essentially equal to the 
most probable log frequency of the synthetic data, with 
a Lomb peak not so uncommon either (corresponding 
approximately to a 10% significance level). 

Aftershock sequence 30 of the Santa Barbara Island 
earthquake is found to exhibit the strongest log-periodic 




t (day) 

Figure 10. Cumulative number of events at the op- 
timal magnitude cutoff 4.5 for the 1933 Long Beach 
earthquake aftershock sequence. The empirical data are 
represented by circles connected by lines to guide the 
eye. 
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Figure 11. Residual oscillations of the fit of Figure 
10 by a simple power law are themselves fitted by a 
simple cosine function (solid line). The empirical data 
are represented by circles connected by lines to guide 
the eye. 
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Figure 12. Spectral analysis of the oscillations shown 
in Figure 11 using the Lomb Periodogram for unevenly 
sampled data. 




Figure 13. Distribution of the log angular frequency ui 
and peak height of 200 synthetic data sets and compar- 
ison to the real sequence identified as the intersection 
of the two lines. 
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Figure 14. Log-periodic oscillations around the cu- 
mulative number of events of the Santa Barbara Island 
aftershock sequence 30. 



oscillations, shown in Figure 14 and 15. However, from 
10 synthetic sequences using the parameters from the 
real sequence, we observed four better log-periodic os- 
cillations at similar log frequencies. We also find that 
the amplitude of oscillations of the real sequence is close 
to that of the synthetic sequences, and this is a general 
property found for all the other sequences. From more 
extensive simulations with 200 synthetic sequences we 
confirm that these oscillations shown in Figure 14 are 
near the most probable region of the distribution of the 
log frequencies in the synthetic time series. 

To summarize these results, let us first discuss the 
time series obtained by using the magnitude cutoffs cho- 
sen by Kisslinger and Jones [1991]. We find that only 
10 sequences among the 27 analyzed sequences have log 
frequencies outside the most probable range of synthetic 
log frequencies. Among them, seven (sequence 11, 13, 
20, 22, 28, 29, 30) have log frequencies slightly higher 
than the most probable synthetic frequencies and might 
thus be considered as potential candidates for a genuine 
physical origin. The three other cases have log frequen- 
cies that are lower than found for most of the synthetic 
tests, a signature of a trend or large-scale aliasing. With 
respect to the statistical significance of the peaks ob- 
tained from the Lomb periodogram, only six sequences 
(10, 18, 23, 27, 31, 36) have spectral peaks higher than 
the most probable synthetic peak heights, and none of 
them overlap with the seven candidates selected above. 

For time series constructed with the best magnitude 
cutoff we find that eight of the log frequencies are higher 
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Figure 15. Normalized Lomb periodogram of the os- 
cillations in Figure 14. The horizontal line marks sig- 
nificance level of 0.01. 



than those found at the magnitude cutoffs of Kisslinger 
and Jones [1991], nine are lower, and ten are equal. 
However, the log frequencies from the corresponding 
synthetic data sets are also shifted accordingly. Five 
sequences (5, 11, 15, 23, 30) have log frequencies higher 
than the most probable range of the frequencies from 
the synthetic data sets, six sequences (1, 13, 14, 22, 27, 
36) have peak higher than those from synthetic data 
sets, but there is no overlap between these two sets 
of five and six sequences, respectively. No significant 
improvement was thus gained by optimizing the mag- 
nitude cutoffs for more regular log-periodic oscillations. 
These results are summarized in Tables 1 and 2. When 
our choice of the best magnitude cutoff is the same as 
that of Kisslinger and Jones, the sequence is excluded 
from Table 2 to avoid repetition. 

We thus find that the 27 aftershock sequences in 
southern California analyzed by the cumulative dis- 
tribution method exhibit clear log-periodic structures, 
which are, however, not strongly statistically signifi- 
cantly different from those generated in synthetic time 
series. Indeed, synthetic series generated using the pure 
Omori's power law exhibit similar structures that re- 
sult from the noise reddening mechanism explained in 
sections 2.1 and 2.2. Optimization over magnitude cut- 
off for more regular oscillations does not significantly 
improve the results. Wo observed that log-periodic os- 
cillations with angular log frequencies are found in the 
interval [0.8,6.2] in most sequences with various regu- 
larity. When compared with synthetic aftershock se- 



quences, these frequencies and the power spectrum at 
these frequencies are found among the most probable 
values in synthetic sequences with the same parameters 
and analyzed in the same way. The amplitudes of the 
oscillations in the real sequences are almost equal to 
those of the synthetic data. 

Another natural quantity to investigate in order to 
qualify Omori's law (1) is to take the numerical deriva- 
tive of the cumulative number M of aftershocks with 
respect to the variable Int. For a pure Omori's law with 
exponent p = 1, this log derivative should be constant. 
Any deviations from a constant provides a measure of 
deviation from Omori's law. We estimate the logarith- 
mic derivative by 



dM 



d\n{t) 



ln(ij+i) - In {ti-i) 



(18) 



Owing to the intermittent nature of the data, a moving 
average has been applied to the derivative. An example 
of the result is shown in Figure 16. Rather than a con- 
stant that would qualify Omori's law, we observe quite 
strong approximately log-periodic signatures. However, 
a test on a synthetic data sequence, using the same 
number of points and time window as well as the same 
Omori's p estimated by Kisslinger and Jones [1991], 
gives a similar signal, as shown in Figure 17, namely 
approximately equidistant peaks of fluctuations. The 
log-periodic signatures are generated from the noise and 
are probably not genuine. In fact, the most obvious dif- 
ference between the true data and the synthetic data 
is that the oscillations in the synthetic data have been 
shifted down to shorter timescales. This is presumably 
due to an imprecise estimation of the value of the ex- 
ponent p. 

2.4. Cumulating and Log-Periodicity 

The essential property of the logarithmic sampling 
is a flat noise spectrum (in the log variable) before in- 
tegration. Integration then leads to the appearance of 
a low- frequency peak, as discussed in Section 2. Dif- 
ferent samplings will give rise to different noise spectra, 
which, when convoluted with the integration, might well 
not create spurious log-periodicity, or even destroy an 
existing one. We now give an example of the latter 
situation. 

We have generated 10 data sets with 30 points each 
using 

y(f) =t-0-5[l + 0.1cos(71nt)] , (19) 

where the sampling was random in the sense that 
the spacing between two consecutive points was cho- 
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Table 1. Log-Periodic Signals in the Oscillations Around the Modified Omori's Power Law of the Cumulative 
Number of Aftershocks Using Kisslinger and Jones' Magnitude Cutoff 



seq 


mm 














pkhto^ 


pkhti^ 


pkht'' 


Pr(Chance)* 


01 


6.3 


1, 


.0740 


2, 


.7077 


1 


,7175 


14.4865 


35.6532 


14 


0.63 


02 


5.1 


1, 


.4913 


4, 


,2834 


2, 


.1367 


1.6963 


6.9105 


5 


0.1974 


03 


6 


1, 


.8340 


2. 


.5803 


2 


.4112 


4.0659 


7.8017 


5.5 


0.4300 


05 


6.9 


2, 


.2869 


2, 


,8769 


1 


.8280 


3.7103 


8.1253 


8 


0.1548 


06 


5.9 


0, 


.7297 


1, 


,7442 





.8516 


8.0163 


17.4935 


9 


0.8495 


09 


6.2 


1, 


.1510 


1. 


,4220 


1 


.2025 


22.4059 


29.2246 


28 


0.3152 


10 


5.5 


1, 


.5570 


2, 


,6391 


1 


.6537 


9.7964 


22.9913 


37 


0.005 


11 


6 


0, 


.8838 


1, 


,8858 


2 


.3547 


13.0106 


25.0258 


16 


0.0606 


12 


5.9 


1, 


.2057 


2, 


,9265 


1 


.5851 


11.9195 


26.1618 


18 


0.6 


13 


5.7 


1 


.2308 


1 


,7628 


2 


.1031 


3.4584 


7.6336 


5.8 


0.2479 


14 


7.7 


1, 


.1584 


1, 


,3052 


1 


.1763 


31.2361 


58.8183 


44 


0.4724 


15 


6.2 


0, 


.8841 


2, 


,1099 


2 


.0406 


9.4233 


22.2794 


16 


0.0414 


17 


5.6 


1, 


.6172 


2, 


,4806 





.9766 


1.4954 


3.7202 


3.7 


0.28 


18 


6.5 


0, 


,9946 


2, 


,3998 


0, 


.9688 


11.3033 


18.0648 


26 


0.0118 


20 


6.6 


0, 


.6763 


1, 


,6025 


2 


.0684 


19.4909 


46.3860 


25 


0.0863 


22 


5.1 


1, 


.1030 


2, 


,2221 


3 


.6476 


38.4192 


47.8631 


37 


0.0125 


23 


5.8 


1, 


.1490 


2, 


,6941 


1 


.4775 


9.1692 


17.2860 


18 


0.1907 


24 


5 


1, 


.1602 


2, 


,2417 


1 


,2623 


13.9662 


28.9842 


17 


0.7714 


25 


5.2 


1, 


.6197 


2, 


,9343 


1 


,9517 


31.0554 


50.2901 


42 


0.3984 


27 


5.5 


1, 


.0248 


2, 


,3511 


1 


,1171 


4.7256 


15.4207 


16 


0.0914 


28 


6.1 


1, 


.9769 


5, 


,0416 


6 


,2079 


2.6342 


9.3955 


8 


0.0102 


29 


5.7 


2, 


.1068 


3, 


,3735 


3 


,4648 


6.8699 


18.2430 


14 


0.0278 


30 


5.3 


1, 


.1045 


2, 


,2538 


2 


,6459 


22.4250 


45.6508 


32 


0.0408 


31 


6.5 


0, 


.7557 


1, 


,8494 


1 


,1278 


39.3589 


68.9044 


75 


0.0811 


32 


5.8 


1, 


.1349 


2, 


,5303 


1 


,8404 


31.0060 


53.6768 


36 


0.3254 


33 


5.9 


1, 


.2886 


2, 


,3381 


1 


.7686 


157.7374 


430.4343 


400 


0.0950 


36 


6 


1, 


.0569 


2. 


,0341 


1 


.4477 


25.5461 


51.6242 


62 


0.0308 



"The number of the data sets used by Kisslinger and Jones [1991]. 
'Main shock magnitude. 

"^The left- wise angular frequency at half the maximum of the distribution of angular frequencies from synthetic aftershocks. 
''The right- wise angular frequency at half the maximum of the distribution of angular frequencies from synthetic after- 
shocks. 

"^Prom the real aftershock sequence. 

^Thc loft- wise angular frequency at half the maximum of the distribution of the peak heights of the most significant Lomb 
peak from synthetic aftershocks. 

'The right- wise angular frequency at half the maximum of the distribution of the peak heights of the most significant 
Lomb peak from synthetic aftershocks. 

''Peak height of the Lomb periodogram for the real aftershock sequence. 

'Probability of the log-periodicity in real aftershocks being chance. 
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Table 2. Log-Periodic Signals in the Oscillations Around the Modified Omori's Power Law of the Cumulative 
Number of Aftershocks Using Other Magnitude Cutoff^ 



seq # 




Wo 








pkhto 


pkhti 


pkht 


Pr(Chance) 


01 


1 


.7518 


3, 


.6436 


2.4358 


3.3245 


8.6912 


9 


0.1181 


03 


1 


.0438 


3, 


.3507 


1.9949 


6.0207 


12.8901 


11 


0.0417 


05 


1 


.0686 


2, 


.8314 


11.9837 


4.5386 


10.1490 


8 


0.0068 


11 


1 


.0552 


2. 


.3074 


2.3547 


9.5031 


17.9494 


16 


0.0606 


13 


1 


.2473 


2, 


.8182 


2.1031 


3.3596 


10.3787 


11 


0.0213 


14 


1 


.0120 


2, 


,1585 


1.04 


10.1416 


22.0273 


28 


0.0357 


15 





.9610 


2, 


.0448 


2.2577 


12.5354 


34.1214 


30 


0.0132 


17 


1 


.1423 


3, 


.4385 


1.0949 


1.3481 


6.9516 


6 


0.375 


18 


1 


.0009 


2, 


.5511 


2.5486 


15.9492 


34.7236 


30 


0.0301 


20 


0, 


.7923 


1, 


.8448 


1.3063 


52.6478 


118.4051 


95 


0.1019 


22 


1, 


.7097 


3, 


.0411 


1.0599 


50.9766 


71.1400 


95 


0.005 


23 


1 


.4364 


2, 


.4938 


2.6761 


19.0478 


56.8038 


40 


0.0455 


24 


1 


.3325 


2, 


.6773 


1.2826 


7.5998 


17.1631 


15 


0.3005 


25 


1 


.6782 


3, 


.6164 


1.7050 


49.2585 


71.3087 


10 


0.005 


27 





.9257 


2, 


.2607 


1.0769 


7.9868 


22.5240 


24 


0.0655 


28 


1 


.7700 


3, 


.8280 


2.4337 


42.6775 


78.4212 


12 


0.005 


29 


1 


.8001 


3, 


.1963 


2.4748 


12.6694 


26.2655 


25 


0.1224 


30 


0, 


.9557 


2, 


.4180 


2.6459 


16.8089 


33.8548 


31 


0.0408 


32 


1 


.1580 


2, 


.7235 


1.6682 


13.4569 


33.1647 


22 


0.3408 


33 





.9972 


2, 


.8735 


2.2487 


3.4352 


9.3050 


9 


0.0588 



"Meanings as in Table 1. 
"Meanings as in Table 1. 
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Figure 16. Aftershock sequence following the Owens 
Valley M5.8 earthquake on October 14, 1978. Deriva- 
tive (denoted "log derivative" ) of the cumulative num- 
ber of aftershocks with respect to the logarithm of time 
as a function of In (t) is shown. 
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Figure 17. Same as Figure 16 for a synthetic data 
set corresponding to the aftershock sequence shown in 
Figure 16. The same number of points and time interval 
as in Figure 16 were used. The p value was p « 1.249 
as listed by Kisslinger and Jones [1991]. 




Figure 18. Typical artificial data set of 30 points gen- 
erated with Eq. (19) with a random sampling such that 

the spacing between two consecutive points was chosen 
from the interval [0, 20] with a uniform probability. 



sen from the interval [0, 20] with uniform probabil- 
ity. From the Lomb periodogram of the log derivative 
[lnt/(ti+i) - lny{ti)]/[lnti+i - Inti] shown in Figure 19, 
we can clearly identify the log-periodic component in 
the signal and its log frequency / = 7/27r~l.l with a 
level of significance better than 0.995 (see Figure 22). 

In contrast, the log-periodic signal completely disap- 
pears in the cumulative distribution as seen in Figure 
23. This not only shows that a strong log-periodic signal 
will be destroyed by "integration" but also that with- 
out noise, only two to three oscillations are needed to 
qualify log-periodicity. 

It is well-known that an integration of sampled data 
corresponds to a low-pass filter. What we have seen 
here is that this low-pass filtering extends into the log 
frequency domain with surprising efliciency. 



3. Mciximum Likelihood Analysis of the 
Aftershock Sequences 

The previous analyses illustrate the danger of ap- 
parently innocuous manipulations of the data, such as 
calculating cumulative distributions of aftershock se- 
quences. This suggests that a method that gives di- 
rect access to the local Omori's exponent may be better 
suited for the luibiased detection of log-periodic struc- 
tures. Indeed, the most direct quantification of depar- 
ture from Omori's law is to look at the logarithmic 
derivative of the total number N of events until time 
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Figure 19. The logarithmic derivative of the data set 
in Figure 18. 
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Figure 20. The cumulative distributions of the data 
set in Figure 18. 
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Figure 22. The Lomb periodogram of the derivative 
in Figure 19 at / = 7/27r rj 1.1. 
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Figure 23. The Lomb periodogram of the derivative 
in Figure 21. 



Figure 21. The logarithmic derivative of the cumula- 
tive distributions in Figure 20. 
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dlnN{t) 
dint 



-Pit) 



which defines the local {i.e., at time t) Omori's exponent 
Pit). 

The maximum likelihood method is nonparametric 
and gives the most probable vahic of p given the data. 
In practice, we perform the maximum hkelihood deter- 
mination in a moving window from t to tjj. Moving t 
and tjj over the time axis allows us to capture local fluc- 
tuations in the Omori's p and provides us with a method 
whereby we can directly estimate any log-periodic os- 
cillations [Lee, 1999]. 

3.1. Maximum Likelihood Estimation of 
Omori's Exponent 

The maximum likelihood method estimates the p 
value in the window from t to tu by maximizing the 
probability that the particular temporal distribution of 
the N events in that window results from Omori's law 
with this particular p value. The most likely p value is 
then given by the implicit equation 



1 t^rT^ \nt - tP-Hiitu 



p-1 



+ 



t'. 



u 



tP- 



(20) 



where the t„ are the time occurrences of the N events 
between t and tu and 



(Int) 



N 



1 ^ 



N 



(21) 



n=l 



An expansion in powers of {p — l)lntu/t gives the ex- 
plicit solution 



12- 



(Ini) 



(ln^r + 1 



with a typical error 



/ 12 
{N-l){\nify 



(22) 



(23) 



obtained from the negative Hessian of the log likelihood. 
Only the exponent p is estimated from the formulas. 
Expression (22) gives an extremely good approxima- 
tion to the exact solution of Eq. (20), since the first 
nonzero correction is proportional to the fourth power 
« 10~^ [{p — l)\ntu /t]'^ of the expansion parameter. 
We have verified that the values p given by (22) are in- 
distinguishable from the exact solution of (20) for the 



range of p values of interest for aftershocks. Since the 
error a is inversely proportional to the inverse of the 
square of the number N of events in the time window, 
probing smaller timescales is bought at a price of larger 
errors. 

3.2. Generation of Log-Periodicity 

Since the precision on the determination of p is essen- 
tially controlled by the number of points in the run- 
ning window, it seems natural to vary tu as a function of 
t such that the number TV of points in the running win- 
dow remains fixed. Any other specification deteriorates 
the homogeneity of the determination of the exponent 
along the time axis. However, this procedure turns out 
to generate strong spurious log-periodicity as well. This 
is illustrated in Figure 24 showing the local p value in a 
running window with N = 100 events determined from 
a synthetic aftershock series generated using (16). The 
spectral analysis of this data set shown in Figure 25, 
confirms the statistical significance of the log-periodic 
signal. 

By construction, the synthetic aftershock sequence 
obeys exactly Omori's law with p= 1, and any observed 
structure decorating this simple power law reflects ran- 
dom noise and finite size effects. We trace back the 
log-periodicity observed in Figure 24 to the constraint 
that the number of points N in the analyzing window 
is fixed. Indeed, for Omori's law with p = 1, this cor- 
responds exactly to a sampling that is uniform in log 
time. As for the cumulative method, the most probable 
component of the noise has a significant spectral power. 

3.2.1. Log frequency as a function of window 

size. To demonstrate the origin of log-periodicity in 
the maximum likelihood method with a fixed number 
of points in the running window, we generate 1000 syn- 
thetic aftershock series with 570 events and p = 1.282, 
corresponding to the values observed for the aftershock 
sequence of Long Beach 1933 documented by Kisslinger 
and Jones [1991]. For each of the 1000 synthetic se- 
quences we have determined p{t) in running windows of 
size TV = 12, 25, 50, and 100 events. For a fixed win- 
dow size A^ we took the Lomb periodogram of each p{t) 
and then averaged them over the 1000 realizations. We 
thus obtain averaged Lomb periodograms for the four 
different values A^ = 12, 25, 50, and 100, shown in the 
Figures 26 and 27. It is clear that the dominating log 
frequency is a function of window size N . 

The approximate lognormal shape of the Lomb peri- 
odograms shown in Figure 27 is remarkably consistent 
across the different values of the window size and is 
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Figure 24. Omori's exponent p as a function of tlie 
logaritlim of time in a running window oi N = 100 
events for a synthetic aftershock time series generated 
using (16) with p = 1.0, extending from time 1 to 1000 
and containing ~ 500 events. Log-periodicity is clearly 
evident to the naked eye. 
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Figure 25. Lomb periodogram of Omori's exponent p 
shown in Figure 24 exhibiting a highly significant peak. 
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Figure 26. Average Lomb periodogram over 1000 
synthetic aftershock realizations with 570 events and 
p = 1.282 of the function p{t) calculated with the max- 
imum likelihood method in running windows of size 
N = 12, 25, 50, and 100 events, respectively. 
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Figure 27. Same as Figure 26 using a log scale for the 
log frequencies illuminating the very similar approxi- 
mate lognormal structure of the Lomb periodograms. 

characterized by approximately the same standard de- 
viation found close to 0.45 for all cases. A plot of the 
logarithm of the log frequency versus the logarithm of 
the window size N is given in Figure 28 showing that 
the log frequency is inversely proportional to N. For a 
point process with Omori's exponent close to 1, keeping 
a constant window size N is approximately equivalent 
to having windows with a constant value of tu jt where 
tu and t are the upper and lower times of the window: 

where N is window size. This reveals that using a fixed 
number of observations in each window corresponds to 
log sampling, a mechanism already identified as a source 
of log-periodicity. 

Of course, for a single realization, fluctuations scram- 
ble the clear-cut result (24) . This is illustrated in Figure 
29 which tests the dependence of the most probable log 
frequencies of the p{€) function as a function of window 
size N varying from 10 to 500. Figure 29 plots the most 
probable log frequency as a function of window size N 
for three cases. The first one corresponds to a synthetic 
sequence of 2000 events with exponent p = 1.0. The 
second one is for a real aftershock sequence, the North 
Palm Springs sequence of 1986, which contains 1437 
events with p k, 1.14. The last example is a numerical 
simulation of 2000 events from the p = 1 stress corrosion 
model of aftershocks studied [Lee and Sornette, 2000] 
with conservation parameter 7 = 0.15, stress corrosion 
exponent a = 15, and stress drop parameter f3 = 0.999, 
for which log-periodicity decorating Omori's law is very 
strong and thus strikingly visible directly on the p{t) 
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Figure 28. Dependence of the most probable log fre- 
quency (in logarithmic scale) corresponding to the max- 
imum of the peaks in Figures 26 and 27 as a function 
of the logarithm of the window size N. The dashed line 
has a slope —1. 
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Figure 29. Log frequency as a function of window 
size N for a synthetic sequence, a real sequence (North 
Palm Springs of 1986 which contains 1437 events), and 
a sequence from the stress corrosion model studied by 
Lee and Sornette [2000]. 
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Figure 30. p{t) as a function of the logarithm of time 
for a synthetic sequence generated with p = 1.5, 2000 
events. 



function (not shown). At large window sizes, multiple 
periods of log oscillation arc included within a window, 
and Omori's exponent p starts to average out giving a 
very low value for the log frequency. At very small win- 
dows, small-scale variations dominate yielding anoma- 
lously high log frequencies. The optimal window sizes 
are small enough to cover less than a whole period of 
oscillation in order to resolve the oscillations but large 
enough not to be dominated by the small-scale varia- 
tions. For the numerical simulation of the stress cor- 
rosion model of aftershocks studied by Lee and Sor- 
nette [2000], we observe a well-defined plateau for the 
log frequency as a function of window size N, extend- 
ing from N = 60 to above 400, qualifying a genuine log- 
periodicity : the theoretical value of the log frequency is 
0.41, [Lee and Sornette, 2000], which is consistent with 
the measured values. In contrast, the log frequency 
dependences of the synthetic Omori's law and of the 
North Palm Springs sequence are very similar to each 
other and show an overall decay as predicted by (24), 
decorated by fluctuations. 

3.2.2. Log frequencies as a function of Omori's 
exponent p. Varying the value of Omori's exponent p 
away from 1 yields quantitatively different results. The 
functions p{t) measured from synthetic time series gen- 
erated with p = 1.5 and p = 0.5 are given in Figures 
30 and 31. Three features can be observed: (1) the 
log frequency has increased in both cases; (2) the os- 
cillations are no longer truly log-periodic since the log 
frequency decreases with time for p = 1.5 and increases 
with time for p = 0.5; and (3) the amplitude of oscilla- 
tions changes with time in both cases. To recover reg- 
ular oscillations, we see from (16) that one would need 
to plot p{t) as a function of t^~P instead of Int, which is 
the correct variable only for series with constant p= 1. 
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Figure 31. p{t) as a function of the logarithm of time 
for a synthetic sequence generated with p = 0.5, 2000 
events. 
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Figure 32. Average log frequency as a function of p 
value for synthetic aftershock sequences. 

Notwithstanding its variation with time, we deter- 
mine the log frequency of the dominant peak by aver- 
aging over the Lomb periodograms of 100 realizations 
of synthetic data. The dependence of the average log 
frequency as a function of the exponent p used in gen- 
erating synthetic aftershock sequences is presented in 
Figure 32. For each p value the average frequency has 
been averaged over 100 realizations. As p differs from 1, 
the sequence exhibits higher and higher log frequencies. 
Figure 32 recovers the results already shown in Figure 
8 for the cumulative method. 

We have attempted to use the insight provided by 
Figures 28, 29, and 32 to distinguish the log-periodic 
signatures observed in real aftershock sequences from 
those obtained in synthetic data sets with the maximum 
likelihood approach. We estimate the log frequency in 
the real aftershock sequences as shown by the arrows 
in Figure 29 on a plateau in a central region for the 
window size which seems to provide a reasonable com- 
promise as explained above. Each of the real aftershock 
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Figure 33. Average log frequency as a function of 

Omori's exponent p in real and synthetic sequences. 
The confidence bounds are seen from the size of the 
fluctuations of the synthetic data. 

sequences is characterized by its average p value de- 
termined by Kisslinger and Jones [1991]. We use this 
p value to generate 1000 synthetic sequences with the 
same number of events. We record the log frequency of 
the largest peak of the Lomb periodogram of the func- 
tion p{t) for each of these 1000 synthetic sequences and 
then average it. Spanning the different real aftershock 
sequences allows us to span a relatively large interval 
of p values. We thus represent in Figure 33 the most 
probable log frequency as a function of average p value 
obtained from the set of real aftershock sequences and 
their corresponding synthetic ensemble. As in Figure 
32, the synthetic log frequency average is increasing for 
p> I (notwithstanding the averaging over 1000 realiza- 
tions, fluctuations remain apparent). In contrast, the 
most probable log frc;quency obtained for the real after- 
shock sequences do not seem to be noticeably dependent 
on the average p. This diff'erence suggests that log- 
periodicity in real aftershock sequences could be real. 
However, the diS'erence is weak and not clearly above 
noise level. By itself, this result is suggestive but not 
enough to substantiate the claim of log-periodicity in 
real data [Lee and Sornette, 2000]. 

3.3. Can Log-Periodicity Be "Saved" for 
Aftershocks? 

3.3.1. Fixed time window size. In view of the 
spurious log-periodicity generated when using a con- 
stant number of points in a moving window which is 

equivalent to log sampling, we redo all calculations us- 
ing the maximum likelihood approach with windows of 
flxed duration (for instance one-day windows). These 
calculations are performed on all sequences documented 
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Figure 34. Comparison of the local p(t) as a function 
of tlic logarithm of time for event 29 (Westmoreland 
aftershock sequence) in the list of Kisslinger and Jones 
[1991] using the maximum likelihood method with a 
running window of fixed time duration of 0.1 day and 
with a fixed number of aftershocks (log sampling). We 
observe structures that are common to the two meth- 
ods: this may be the signature of a real log-periodicity. 



by [Kisslinger and Jones, 1991] as well as on the more 
recent Loma Prieta, Landers, and Northridge after- 
shock sequences. As typical examples, we show the 
1981 Westmoreland and 1989 Loma Prieta aftershock 
sequences. We choose the Westmoreland as an exam- 
ple with p > lA and also with fluctuations greater than 
one standard deviation {2a) from the mean equal to one 
(Figure 34). Loma Prieta is chosen as an example of a 
more recent sequence with a more complete record. It 
also exhibits fluctuations in the local p{t) greater than 
2a. Using a running time window which is fixed to 0.1 
days, the Omori's exponent p for the Westmoreland and 
for Loma Prieta aftershock sequences appear to fluctu- 
ate log-periodically with a log frequency around 1. The 
solid lines are for a constant number of events in the 
running window, while the dashed lines are for fixed 
time windows. In the two examples shown here and 
for a majority of the real aftershock sequences studied 
here, peaks are found at log frequencies near 1 (see the 
peaks in solid and dashed lines close to the abscissa 1 
in Figure 35). 

Table 3 lists the top three local peak max;ima in fre- 
quency with their associated spectral power from the 
Lomb periodogram for all analyzable real aftershock se- 
quences, using the maximum likelihood method with 
fixed time window. The sequence numbers are as given 
by Kisslinger and Jones [1991]. The amplitude of the 
fluctuations in p values in running time windows as a 
function of time is given in integer units of the stan- 
dard deviation around the average calculated for each 
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Figure 35. Comparison of the Lomb periodograms 
corresponding to the two signals shown in Figure 34. 
We observe a common peak at a log frequency close to 
1.2. 
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Figure 36. Same as Figure 34 for the Loma Prieta 
aftershock sequence. 
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Figure 37. Same as Figure 35 for the Loma Prieta 

aftershock sequence. We also observe a common peak 
at a log frequency close to 1 and a slightly smaller one 
at approximately half log frequency. 
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Table 3. List of the Top Three Local Peak Maxima in Frequency With Their Associated Spectral Power From the 
Lomb Periodogram for All Analyzable Real Aftershock Sequences, Using the Maximum Likelihood Method With 
Fixed Time Window" 



Sequence Amplitude(a) 


window (days) 


/l 


/2 


h 


poweri 


powers 


powers 


1 


2 


150 


4, 


.18 


3 


.61 


2 


.60 


13.74 


11.84 


8.42 


5 


1 


1 





.71 


1, 


.18 


1 


.65 


4.68 


2.52 


2.14 


6 




1 


1 


.61 


1, 


.80 


1 


.38 


2.65 


2.58 


2.26 


7 


1 


100 


1 


.90 


0, 


.79 


1, 


.57 


5.67 


5.16 


5.11 


8 


1 


10 





.56 


1, 


.48 


4, 


.20 


5.85 


3.42 


3.34 


9 


1 


1 


4 


.01 


4, 


.95 


1 


.50 


2.83 


2.63 


2.61 


10 




1 


1 


.79 


1, 


.47 


4, 


.20 


3.13 


2.16 


1.75 


11 




1 


2 


.35 


2, 


.63 





.76 


2.00 


1.69 


1.53 


12 




1 





.63 


3 


.14 


1 


.84 


4.27 


2.23 


1.48 


13 


1 


1 


1, 


.48 


1, 


.55 


2 


.48 


3.01 


2.95 


2.24 


14 


2 


100 


1 


.59 


3, 


.65 


2, 


.06 


16.73 


13.43 


5.49 


15 




1 





.63 


0, 


.89 


2 


.99 


3.37 


2.64 


2.54 


18 




1 


0, 


.52 


4, 


.39 


2 


.87 


4.57 


3.82 


3.78 


19 




10 


0, 


.66 


1, 


.33 


2, 


.78 


3.92 


1.38 


1.31 


20 




1 





.96 


2, 


.82 


2 


.91 


3.91 


2.94 


2.86 


21 




10 


0, 


.61 


1, 


.22 





.84 


2.18 


2.11 


1.60 


22 


1 


1 


2 


.47 


2, 


.30 





.69 


4.79 


4.69 


4.64 


23 




1 


1 


.19 


0, 


.69 


1 


.88 


4.70 


1.38 


1.17 


24 


1 


1 


1, 


.31 


1, 


.07 


1 


.45 


4.96 


3.99 


3.66 


25 


1 


1 
i 


2 


.53 


2, 


.22 


1 


.99 


^^ /I n 
D.4y 




4.U0 


26 


1 


1 





.70 


4, 


.13 


3 


.03 


22.82 


18.04 


11.46 


27 




1 


0, 


.76 


4, 


.34 





.56 


2.56 


2.16 


2.01 


28 


2 


1 


1 


.25 


1, 


.47 


2, 


.36 


7.45 


7.07 


5.86 


29 


1 


0.1 


1 


.30 


3 


.84 


4, 


.87 


13.07 


6.57 


6.08 


30 




1 


0, 


.66 


3 


.55 





.96 


5.61 


4.69 


3.07 


31 


2 


1 


4, 


.06 


3, 


.45 


1 


.72 


15.44 


5.87 


5.37 


32 




1 





.53 


1, 


.38 


1 


.68 


6.62 


5.12 


4.87 


33 


1 


1 


1 


.41 


3 


.57 


2 


.70 


30.14 


23.42 


16.79 


34 


2 


10 





.58 


3, 


.77 


1 


.50 


19.68 


15.43 


7.41 


35 




1 





.55 


4, 


.85 


1 


.65 


3.94 


3.08 


2.31 


36 




1 


1, 


.42 


1, 


.65 


2 


.22 


8.89 


8.42 


7.70 


Loma Prieta 


2 


0.1 





.86 


4, 


.83 


4, 


.55 


15.91 


12.59 


9.53 


Landers 


2 


1 


4, 


.06 


4. 


.85 





.97 


19.14 


17.11 


16.09 


Nort bridge 


2 


10 





.60 


1. 


.04 


4, 


.20 


50.54 


49.38 


19.36 



"The sequence numbers are as Kisslinger and Jones [1991]. The amplitude of the fluctuations in p values in running time 
windows as a function of time is given in integer units of the standard deviation around the average calculated for each 
sequence. 

"The sequence numbers are as Kisslinger and Jones [1991]. The amplitude of the fluctuations in p values in running time 
windows as a function of time is given in integer units of the standard deviation around the average calculated for each 
sequence. 
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sequence. 

3.3.2. Lomb averaging. The case for log-periodi- 
city in earthquake aftershocks is quite weak. We have 
shown that the maximum likehhood method with run- 
ning windows of a fixed number of events generates 
spurious log-periodicity with a log frequency inversely 
proportional to the window size. Furthermore, the mea- 
sured frequencies do not significantly deviate from those 
found from 1000 synthetic sequences. We moved to the 
maximum likelihood method with a fixed time window 
to remove the log sampling effect. The results were 
slightly more encouraging. We have also tried to an- 
alyze fluctuations of the difference between the cumu- 
lative number of aftershocks and its smoothed Omori's 
law. Although the structures generated by this method 
were in good agreement with the p value structure from 
the maximum likelihood method, the fact that both 
methods give similarly correlated structures on syn- 
thetic data does not enhance the evidence. 

We now turn to a last investigation using averaging. 
As shown by Johansen et al. [2000b], averaging over 
Lomb periodograms provides in principle a powerful 
method to retrieve the presence of log-periodicity if the 
log frequencies are similar from sequence to sequence. 
This is due to the fact that the Lomb periodograms 
are not sensitive to the phase of the log-periodic oscil- 
lations, and we do not run the risk of scrambling the 
log-periodic oscillations due to their random phases, 
as woiild occur during a direct averaging of the p{t) 
functions for instance. See also Johansen and Sornette 
[1998] for a general discussion of the caveats of ensemble 
averaging. 

Figure 38 shows the average Lomb periodogram ob- 
tained by using all real aftershock sequences. It should 
be compared to Figure 39 obtained by averaging the 
Lomb periodograms of 1000 synthetic realizations for 
each real aftershock sequence. We see that the peak is 
differently positioned for the real aftershock sequences 
compared to that of the synthetic. 

Encouraged by this, we have singled out all real af- 
tershock sequences where the fluctuations around the 
average Omori's law of p = 1 was greater than 1 and 2 
a. The result is shown in Figures 40 and 41. We see 
that a stronger log-periodic signal appears in the spec- 
tra of the real aftershock sequences compared to that of 
the synthetic. By this, we observe that the peak culmi- 
nating at 1 of the solid line (real data) is significantly 
sharper than that for the synthetic ensemble. This in- 
dicates that genuine log-periodic structures might be 
most easily identified in sequences with "anomalous" 
high or low p value and large fiuctuations. 



0.6 




2 4 6 8 10 12 14 



Frequency 

Figure 38. Average Lomb periodogram over all real 

aftershock sequences using the maximum likelihood 
method with running windows of a constant number 
of events (log sampling). 
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Figure 39. Same as Figure 38 but for 1000 synthetic 
sequences for each real aftershock sequence with the 
same parameters. 
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Figure 40. The solid line is the average Lomb peri- 
odogram over the subset of 17 real aftershock sequences 
for which the p value fluctuations are greater than la. 
The dashed line is the average over 17,000 synthetic 
aftershock sequences, 1000 for each of the 17 real ones 
with p value fluctuations greater than la. 
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Figure 41. The solid line is the average Lomb pe- 

riodogram over the STibsct of seven real aftershock se- 
quences for which the p value fluctuations are greater 
than 2c7. The dashed line is the average over 7000 syn- 
thetic aftershock sequences, 1000 for each of the seven 
real ones with p value fluctuations greater than 2a. 

4. Discussion 

We have presented a detailed study of log-periodicity 
in relatively sparse data with nonuniform sampling. We 
have focused our attention on the effect of this nonuni- 
form sampling and its interplay with integration in the 
presence of noise, which may lead to synthetically gen- 
erated log-periodicity that may compete with genuine 
signals. Our main result is the discovery that ap- 
proximately logarithmic sampling and integration (per- 
formed either by constructing cumulative functions or 
using a maximum likelihood method) generate strong 
log-periodic structures. This result calls for a reassess- 
ment of previous reports of log-periodicity in rupture 
[Anifrani et al., 1995; Sahimi and Arbabi, 1996; Huang 
et al, 1997; Johansen and Sornette, 1998] and in earth- 
quakes [Sornette and Sammis, 1995; Saleur and Sor- 
nette, 1996; Saleur et al., 1996a, b; Johansen et al., 
1996; Varnes and Bufe, 1996; Lee, 1999] in the spirit of 
the reanalysis [Johansen et al., 2000a] of the evidence 
for criticality and especially log-periodicity in the pre- 
viously reported chemical anomalies that preceded the 
Kobe earthquake. In this recent work [Johansen et al., 
2000a], the ion (CI", K+, Mg++, NO3 and SO J") con- 
centrations of groundwater issued from deep wells lo- 
cated near the epicenter of the 1995 Kobe earthquake 
are taken as proxies for the cumulative damage pre- 
ceding the earthquake. Using both a parametric and 
nonparametric analysis following the guidelines given in 
this paper, the five data sets are compared extensively 
to synthetic time series. The null hypothesis that the 
patterns documented on these times series result from 



noise decorating a simple power law is rejected with a 
very high confidence level. For the other cases of re- 
ported log-periodicity, we plan to present a reanalysis 
in forthcoming publications. 

We have also analyzed in detail the 27 best after- 
shock sequences studied by Kisslinger and Jones [1991] 
and the Loma Prieta, Landers, and Northridge after- 
shock sequences for the presence of log-periodic correc- 
tions to Omori's law. Both the cumulative distribution 
and maximum likelihood methods lead to synthetic log- 
periodicity resulting from the interplay between noise 
reddening and power law scaling. However, it is im- 
portant to distinguish between the two processes. The 
low-pass filtering resulting from the integration used to 
generate the cumulative distribution is very different 
from the averaging process used in the maximum likeli- 
hood method. In the latter, the upper frequency cutoff 
is completely determined by the data window used in 
the estimation process. Most of the structures observed 
on the 27 aftershock sequences can be explained solely 
on the basis of the spurious log-periodicity generated 
by the nonuniform sampling and the integration mech- 
anism. There arc some residual structures in aftershock 
sequences of possible log-periodic origin which cannot 
be fully explained by the synthetic scenario. Hence the 
possibility of a signal resulting from a physically based 
log-periodicity as claimed by Lee [1999] cannot be ruled 
out. Better quality data with an order of magnitude 
more events are needed to resolve this issue. 

A tantalizing question is whether the spurious log- 
periodic oscillations that might appear in the analysis 
of a foreshock sequence can still be used to improve the 
prediction of the main event. This would be equiva- 
lent to "fitting the most probable form of noise," an 
attractive though maybe not workable idea. We hope 
to report on this point soon. 
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